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A  THIRD-ORDER  ANALYTICAL  SOLUTION 
FOR  RELATIVE  MOTION  WITH  A  CIRCULAR 

REFERENCE  ORBIT 


David  L.  Richardson*  and  Jason  W.  Mitchell^ 


Abstract 

From  a  Lagrangian  approach  for  the  development  of  the  relative  motion 
equations  for  a  nonlinear  Hill’s  problem,  it  is  shown  that  the  influence  of  a 
spherical  primary  mass  takes  the  form  of  a  third-body-like  disturbing  func¬ 
tion  expressed  relative  to  the  origin  of  Hill’s  rotating  frame.  The  resulting 
Lagrangian  and  corresponding  equations  of  motion  are  compact  and  can 
provide  an  easily  obtainable  representation  of  the  nonlinear  contributions 
to  the  motion  to  an  arbitrary  order  using  recursion  relations.  The  relative- 
motion  equations  are  expanded  through  third-order  in  the  local  Hill’s  coor¬ 
dinates  and  a  correspondingly  accurate  successive  approximations  solution 
is  developed  to  describe  nonlinear  periodic  motions  in  the  Hill’s  frame. 


INTRODUCTION 

Considerable  attention  has  been  given  to  the  dynamics  and  control  of  formation  flying  satel¬ 
lites  [1, 10,  6].  In  particular,  much  effort  has  been  given  to  determining  the  relative  effects  of 
perturbations  on  satellite  formations.  Frequently,  preliminary  mission  analyses  begin  with 
what  are  commonly  known  as  Hill’s  [4]  or  Clohessy-Wiltshire  [3]  equations.  These  are  the 
linearized  equations  describing  the  relative  motion  of  two  bodies  under  the  gravitational 
influence  of  a  point-mass  central  body.  Traditionally,  these  equations  have  been  used  to 
describe  rendezvous  maneuvers  and  typically  prove  useful  only  for  a  few  orbital  revolutions. 
More  in-depth  analyses  have  included  the  effects  of  central  body  oblateness,  solar  and  elec¬ 
tromagnetic  perturbations,  drag  [11];  and  eccentricity  effects  [5].  Consequently,  it  seems  of 
interest  to  develop  an  analytical  solution  describing  periodic  motions  in  Hill  s  frame  con¬ 
taining  the  most  significant  nonlinear  contributions,  which  influence  the  relative  motion  as 
much  as  the  effects  of  more  well-known  perturbations. 

In  this  paper,  we  describe  a  third-order  analytical  perturbation  solution  developed  by 
successive  approximations  that  provides  a  very  accurate  periodic  solution  to  a  set  of  non¬ 
linear  Hill’s  equations. 

*  Professor,  Aerospace  Engineering  and  Engineering  Mechanics,  University  of  Cincinnati,  P.O.  Box  210070, 
Cincinnati  OH  45221. 

^Visiting  Scientist,  Air  Vehicles  Directorate,  Control  Theory  Optimization  Branch  (VACA),  2210  Eighth 
St.,  Wright-Patterson  AFB,  OH  45433. 
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THIRD-ORDER  NONLINEAR  HILL’S  EQUATIONS 


The  drag-free,  nonlinear  Hill’s  equations  describing  the  motion  of  a  slave  relative  to  a  master 
in  a  circular  orbit  of  radius  R  and  mean  motion  n  are  obtained  from  the  relative-motion 
Lagrangian 

L  =  h-r+  —  +  n^RX,  (1) 

2  p 

where 

r-Xi  +  rj  +  Zk, 

and 

p^  =  P?  +  r‘^  +  2r-R.. 

X,  Y,  Z  form  a  local  dextral  system  of  coordinates  with  origin  at  the  master  as  shown  in 
Figure  1.  The  X  axis  (i)  lies  in  the  direction  of  the  radial  line  outward  from  the  central 
body  (Earth),  Y  is  oriented  in  the  along-track  (velocity)  direction  (j)  and  Z  lies  along  the 
orbit  normal  (k)  in  the  direction  of  the  angular  momentum. 

It  is  algebraically  convenient  to  introduce  nondimensional  coordinates  and  a  nondimen- 
sional  independent  variable  r  via  the  equations 


X  =  Rx,  Y  =  Ry,  Z  =  Rz, 


dt  ”dr‘ 


(2) 


The  new  Lagrangian  for  the  motion  is  the  scaled  Lagrangian  L  — *  £.  =  L/{Ti^R‘^)  trans¬ 
formed  from  Eq.  (1)  along  with  the  substitution,  GE  =  rv^R? ,  obtained  from  Kepler  s  Third 
Law.  The  result  is 

+  X.  (3) 


The  square  root  contribution  is  recognized  as  the  generating  function  for  the  Legen¬ 
dre  polynomials  Pn{~x/Q)  —  Pn{^ / q)  having  argument  —x/g  where  g  is  the  scaled 

distance  from  the  master  to  the  slave  with 


g  —  zk. 

This  produces  the  infinite  series  expansion 

jC  =  y-g'  +  '£i-erPni^/^)- 

^  n>2 


(4) 


2 


Using  the  above  Lagrangian,  the  equations  of  motion  for  the  slave  relative  to  the  master 
are  obtained  through  third-order  in  the  nondimensional  coordinates.  We  find 


x"  -  2j/'  -  3a;  =  -^  +  2a;(2a;^  -  3y^  -  3z^), 

y"  +  2x'  =  3xy  -  ^y  (Ax"^  -y^  -  , 

z"  +  z  =  3xz  -  ^z  (Ax"^  -y^  -  2^)  • 


(5) 


Linear  Solution  and  Constraints 

If  the  nonlinear  contributions  on  the  right-hand  sides  of  Elqs.  (4)  are  ignored,  we  have  the 
well-known  set  of  linear  Hill’s  equations  that  are  often  used  as  the  basis  for  rendezvous  and 
formation  fiying  analyses: 


x"  —  2y'  —  3a:  =  0, 

y"  +  2x'  =  0,  ; 

z"  +  z=^  0. 

The  characteristic  equation  for  the  x^y  system  is 

A2  (a2  +  l)  =  0. 


(6) 


The  repeated  roots,  A  =  {0,0},  give  rise  to  a  linear  solution  that  contains  a  secular  drift. 
We  write  the  general  solution  as 


x{t)  =  ai  -f  as  cos  r  +  a4  sin  r, 

3 

y(r)  =  02  “  o  +  204  cos  r  —  203  sin  r, 
z[r)  =  05  cos  r  +  06  sin  r. 


(7) 


The  integration  constants  Oj  in  the  general  solution  are  related  to  the  components  of  the 
initial  state  vector  {xo,yo,zo,Xo,yo,Zo)'^  through  the  relations 


ai  =  4a;o  -I-  2yo, 

02  =  yo  -  2xf,, 

03  =  -3x0  -  2yo,  I  (3) 

<14  =  Xq, 

O5  =  Zq, 

06  =  z'q. 

We  impose  the  following  constraints  on  the  initial  conditions  to  insure  that  the  linear  solu¬ 
tion  is  periodic  about  the  origin  of  coordinates: 

oi  =  4xo  +  2t/o  =  0, 

02  =  yo  -  2xo  =  0. 
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(9) 

(10) 


The  algebra  encountered  in  our  nonlinear  developments  is  simplified  somewhat  if  we  rewrite 
the  periodic  solution  to  Eqs.  (6)  in  phase/amplitude  form.  To  this  end,  we  have 


x(r)  =  -Acos(t  + 

y(r)  =  2Asm(T  +  (p),  > 

z(t)  =  Bsm(T  +  Ip). 


(11) 


NONLINEAR  APPROXIMATION 


Since  the  periodic  Hill’s  problem  is  the  relative  difference  in  motion  of  two  gravitationally 
non-interacting  two-body  motions  with  identical  semimajor  axes,  we  do  not  expect  the 
frequencies  of  our  nonlinear  periodic  solution  to  be  a  fnnction  of  the  amplitudes  of  the 
motion.  Nor  do  we  expect  to  have  to  impose  any  equality  constraints  on  the  amplitudes 
of  motion  as  was  required  in  a  similar  successive  approximations  derivation  for  periodic 
motions  about  the  collinear  points  in  the  circular-restricted  problem  [8,  9].  Nevertheless, 
to  provide  additional  algebraic  flexibility,  we  introduce  a  frequency  parameter  u>  and  once 
again  change  the  independent  variable  from  r  to  s  via  the  transformation 


A  -  A 

dr  ^ds 


(12) 


The  successive  approximations  solution  to  the  third-order  equations  of  motion  given  in 
Eqs.  (5)  are  the  asymptotic  series 


x(s)  =  Xl  +  X2  +  X3  +  0(4), 
y{s)  =  yi  +  y2  Ays  +  0(4), 
z(s)  =  Zi  +  Z2  +  Z3  +  0(4). 


(13) 


The  system  frequency  u)  is  assumed  to  have  an  asymptotic  development  of 

tjj  =  \  -{■  ui\  4-  (jJ2  A  0(3).  (f4) 


The  subscripted  functions  above  are  taken  to  obey  the  order  of  magnitude  relations 

Xi,yi  =  0(A^),  Zi  =  0(B%  uJi  =  0(A%  (15) 

To  fourth-order,  the  equations  of  motion  are  rewritten  and  separated  into  the  equations 
of  the  first-,  second-  and  third-order.  For  the  first-order,  we  have 

x"  -  2y[  -  3xi  =  0, 
y'{  -t-  2x'i  =  0, 
z'l  -t-  zi  =  0. 

In  these  equations  and  those  following,  primes  denote 

Xl  =  -Acos(s  +  (p), 
yi  =  2i4sin(s  -I-  <p), 
zi  —  Bsm(s  +  Ip). 


I  (16) 

d/ds.  Prom  Eqs.  (11),  the  solution  is 

(17) 
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The  equations  for  the  second-  and  third-order  are 

^2  “  ^2^2  ““  ^^2  “  — 2u;i  {xi  ~  J/i)  ”  2  “  2/i  ? 

2/2  +  2x^2  =  -2^1  (^/i  +  ^i)  +  3xiyi, 

2^2  +  2:2  =  -2a;i;2Ji  +  3a;iZi, 

and 

x'i  -  2?/3  -  3X3  =  -2wi  (x^'  -  -  2u2  {x’(  -  yi)  -  wfx'/ 

—3  (2xiX2  —  yij/2  “  ^1-2^2)  +  2xi  (2xj  —  3j/f  —  3^^) , 

y'^  +  2x'3  =  -2cJi  (y^'  +  x'j)  -  2a;2  «  +  xQ  -  wfyi' 

+3  (xiy2  +  X2yi)  -  2^1  (4x?  -  y?  -  •^i)  , 

Z3  +  Z3  =  -2uiZ2  -  {2u>2  +  OJi)  z'l  +  3  (xi22  +  X2Zx) 

(4x1  -  yf  -  4) . 

respectively, 

Second-Order  Corrections 

To  eliminate  secular  terms  in  the  particular  solutions  to  Eqs.  (18),  the  (expected)  assignment 

uji  —  0  (20) 

is  required.  After  substituting  from  Eqs.  (17),  the  equations  for  the  second-order  corrections 
are 

X2  -  2y2  -  3x2  =  I  (2A^  +  -B^)  -  cos  2u  -  cos  2v, 

y2  +  2x2  =  —  3A^sin2u,  ^  (31) 

3  1 

22  +  22  =  --AB^sm{u  +  u)  +  sin(t;  -  M)j , 

where  the  abbreviations  in  the  trigonometric  arguments  are 

u  =  s  +  <!>,  V  =  s  +  ip. 

Constructing  particular  solutions  is  algebraically  straight  forward.  The  results  are 
X2  =  -^  (2A‘^  +  3"^^  +  i>l^cos2«+  ^B^cos2i;, 

y2  =  -A^sin2u  -  7B^sin2u,  ^  (33) 

4  4 

22  =  ^>lB|sin(w  +  u)  —  3sin(w  —  u)| . 
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Third-Order  Corrections 


Substituting  the  lower-order  solutions  into  the  right-hand  sides  of  Eqs.  (19)  produces 


X3  -  2?/3  -  3X3 


y'i  +  2x3 


4  +  •2=3 


^ cos(u  -  2u)  -I-  A  cos  u 

3  25 

—  -AB^  cos(u  -I-  2v)  — —  A^  cos  3u, 

4  4 

-|aB^  sin(u  -  2t;)  -I-  A  +  2u'2^  sinu 

3  39  q 

-I-  -AB^  sin(tt  -I-  2v)  -  -5- A^  sin  3u, 

8  o 

2Buj2  sinu  -  3A^B  sin(2u  -I-  v). 


> 


(23) 


It  is  seen  that  the  only  choice  for  removing  the  sinu  term  in  the  23  equation  is  to  impose 
the  condition 


UJ2  =  0. 


(24) 


Accordingly,  the  frequency  corrections  are  all  zero  to  this  order  of  analysis,  as  expected. 
Integrating  to  find  a  particular  solution  for  23  produces 


Z3  =  ^A^B  sin(2u  -I-  v). 
8 


(25) 


It  would  seem  that  there  are  various  terms  in  the  coupled  X3,  y'l  equations  that  could 
contribute  secular  results  to  the  solution.  However,  this  is  not  the  case.  If  the  y'l  equation 
is  integrated  once  and  the  result  substituted  into  the  X3  equation,  we  have 


x'l  4-  X3  =  -  AB^  cos(u  -t-  2'u)  -  3A^  cos  3u, 


with  the  particular  solution 


1  3 

X3  =  -AB^  cos(u  -I-  2v)  -t-  -  A^  cos  3u. 

8  o 


(26) 


This  result  allows  us  to  write  an  expression  for  2/3  which  is  integrated  to  give 


j/3  =  _  1^52  2^)  -I-  ^ A^  sin  3u  -1-  ^ AB2  sin(«  -  2v)  -  ^ A®  sin  u.  (27) 
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Figure  2:  ar-direction  comparison  of  pertur¬ 
bation  solution  to  integration  of  third-order 
differential  equations. 


Figure  3:  y-direction  comparison  of  pertur¬ 
bation  solution  to  integration  of  third-order 
differential  equations. 


COMPLETE  SOLUTION  AND  NUMERICAL  TESTING 


We  collate  the  above  results  and  write  our  third-order  successive  approximations  solution 
to  Eqs.  (5)  as  follows: 


x{t)  =  -  a  cos  u  -  ^  (2A‘^  +  2u  +  ^B‘^  cos  2v 

1  3 

+  -AB'^  cos{u  +  2v)  +  -A^  cos  3u, 

So 

111  7 

y^r)  =  2 A  sin  n  +  sin  2u  —  -  J5^  sin  2v  —  sin(w  -h  2v)  + 

3  9 

-h  -AB'^  sin(u  -  2v)  -  -A^  sin  w, 

8  o 

1  3 

z(r)  Bsinv  +  -AR[^sin(w  +  ?;)  —  3sin(t;  —  ^x)j  +  -A^Bsin(2u  A  v), 


}  (28) 


where 


u  =  rA(t>,  V  —  r 


because  a;  =  1,  thus  making  s  =  r. 

Having  developed  a  perturbation  solution,  we  compared  this  to  a  numerical  integration 
of  the  third-order  nonlinear  differential  equations  given  by  Eqs.  (5).  The  fixed  stepsize 
numerical  integration  was  performed  using  a  Runge-Kutta  4/5  method  [2]  spanning  ap¬ 
proximately  one  day  at  an  orbital  altitude  of  500  km. 

The  initial  conditions  were  selected  by  specifying  the  unsealed  amplitudes  and  phases 
used  in  the  perturbation  solution.  We  chose 

3-  20  km,  B*  -  4  km,  ^  -  0,  ^  =  7r/2.  (29) 

The  appropriately  scaled  amplitudes  A,  B  and  (/>,  *0  where  then  used  with  Eqs.  (28)  along 
with  their  derivatives  to  compute  the  test  set  of  initial  conditions. 


7 


Figure  4:  2:-direction  comparison  of  pertur¬ 
bation  solution  to  integration  of  third-order 
differential  equations. 


diff 

amplitude 

avg 

X 

4.1052  X 

+4.3961  X  10-°° 

y 

2.2746  X  10-°® 

-4.7845  X  10-°° 

z 

1.0596  X  10-°® 

-3.1647  X  10-°® 

Table  1:  Perturbation  to  linear  solution  comparison. 


The  position  difference  between  the  third-order  perturbation  solution  and  the  numeri¬ 
cal  integration  can  be  seen  in  Figures  2—4.  Overall,  the  perturbation  solution  appears  to 
represent  a  periodic  solution  to  the  actual  third-order  equations.  The  accuracy  remained 
consistent  with  the  truncation  errors  in  the  successive  approximations  solution.  There  is 
a  small  along-track  secular  runoff  seen  in  Figure  3.  The  likely  source  for  this  behavior  is 
attributed  to  a  small  difference  between  the  semimajor  axes  of  the  slave  and  the  master. 
This  discrepancy  creates  a  small  difference  in  orbit  periods  which  over  a  one  day  simulation 
interval  is  seen  as  a  small  along-track  separation.  For  orbital  configurations  of  similar  en¬ 
ergies,  this  drift  will  remain  small  on  the  order  of  the  truncation  errors  in  the  development 
of  the  equations  of  motion.  All  this  is  reinforced  in  Figure  3. 

For  comparison.  Table  1  shows  the  difference  between  the  linear  solution  and  the  per¬ 
turbation  solution.  Since  the  differences  are  harmonic,  the  amplitude  and  average  of  the 
differences  between  the  linear  and  perturbation  solution  are  provided.  The  amplitude  of 
the  oscillation  indicates  that  the  differences  are  small,  but  not  insignificant  as  a  perturbing 
effect. 

CONCLUSIONS 

In  this  paper,  a  third-order  successive  approximations  solution  describing  relative  motion 
in  the  neighborhood  of  a  circular  reference  orbit  was  presented.  The  equations  of  mo- 
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tion  were  developed  from  a  compact  relative-motion  Lagrangian  for  two  gravitationally 
non-interacting  satellites.  The  resulting  analytical  solution  compared  well  to  a  numerical 
integration  of  the  third-order  nonlinear  differential  equations  over  a  simulation  period  of  one 
day.  The  perturbation  solution  errors  in  the  radial,  along-track,  and  orbit  normal  directions 
remained  on  the  order  of  millimeters,  centimeters,  and  tenths  of  millimeters,  respectively, 
throughout.  These  differences  remained  small  even  in  the  presence  of  slight  differences  in 
orbital  energies.  Discrepancies  were  consistent  with  local  truncation  errors. 
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